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Abstract 

We investigate globally coupled stochastic three-state oscillators, which we consider 
as general models of stochastic excitable systems. We compare two situations:in the 
first case the transitions between the three states of each unit 1 — > 2 — > 3 — > 1 
are determined by Poissonian waiting time distributions. In the second case only 
transition 1 — > 2 is Poissonian whereas the others are deterministic with a fixed 
delay. When coupled the second system shows coherent oscillations whereas the 
first remains in a stable stationary state. We show that the coherent oscillations 
are due to a Hopf-bifurcation in the dynamics of the occupation probabilities of the 
discrete states and discuss the bifurcation diagram. 



1 Introduction 

Coupled dynamic units are of great importance to understand cooperative 
behavior in deterministic dynamical systems. Many investigations on different 
levels of description like reaction diffusion equations, discrete deterministic 
maps, globally connected (mean field) or nonlocal models etc. give nowadays 
a good understanding of temporal and spatial pattern formation in different 
fields of science [1]. 

In the case of a coupled stochastic dynamics, stochastic source terms have to be 
incorporated in the deterministic dynamics. From simulations and equations 
for the moments it was shown that noise may be the origin of nontrivial 
behavior leading to noise induced order at non- vanishing noise intensities [2]. 
In particular, a new type of non-equilibrium phase transitions was found and 
it was shown that the transmission of signals even benefits from the presence 
of noise [3]. 
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Principal insight of noise induced behavior has been achieved by introduc- 
ing stochastic coupled units with few discrete states. For instance, stochastic 
resonance became a useful application in many different contexts after the 
formulation of two-state stochastic systems with periodically driven rates [4] . 
Also coupled chains have explained front motion in fluctuating media and ar- 
ray enhanced stochastic resonance [5]. Noise induced fluxes in sawtooth like 
potentials and media have also been considered with such stochastic systems 
incorporating a small number of possible states [6] . 

Most of the present simple discrete systems consider Poissonian transitions 
between the n states of a unit. But such approximation, yielding a Markovian 
description, require often the introduction of a large number of variables or 
states. Whereas random walks are often considered with different waiting time 
distributions [7] investigations of stochastic discrete dynamical units with, for 
example, delayed feedback, has started only recently [8,9]. 

In [10] we presented the occurrence of global oscillations in coupled three state 
oscillators with fixed waiting time [9] as a prototypical model for coupled ex- 
citable units. Here we compare this model with its Poissonian counterpart. 
We start from a common description for stochastic three state dynamics with 
general waiting time distributions. Whereas the resulting local dynamics char- 
acterized by its stationary density and power spectrum behaves similarly in 
both cases, the coupled cases exhibit a different behavior. In the globally 
coupled network with Poissonian transitions only stationary states can be 
encountered due to the high randomness of the transitions. Contrary, the de- 
layed network shows coherent oscillations [11] and synchronization is found, 
as observed before in numerical investigations of excitable systems [12]. 



2 Single three state units 



Our prototypical model of an excitable system is shown in Fig. 1. The unit 
aims to mimic a single stochastic excitable system. The transition from state 
1 — > 2 is the excitation step and occurs with rate 7. The transitions from 
2^3 and from 3 — > 1 model the return and take place with arbitrary waiting 
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Fig. 1. A single stochastic three-state unit. 
The process 1 ^ 2 is Poissonian, while 
the transitions 2^3 and 3^1 are are 
governed by the waiting time distribution 
1^2 (r) and 1^3 (r). 
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time distributions W2{t) and w^^r), respectively. As output of the unit we 
assign f{t) = 1 if the unit is in state 3 and f{t) — 0, otherwise. 



First we derive equations which govern the dynamics of the probabihties Pi{t), 
i = 1,2,3. that the system is in state i at time t. To this end we introduce 
the probabihty densities Pi{t,T) that the system is in state i at time t and 
has entered it at time t — t. The relation between these probability densities 
and the probabilities Piit) is given by Piit) = J^dr Pi{t,T). The flow Jj(t) 
of probability between state i and i + 1 can be expressed in terms of the 
probability densities Pi(t,T) and the waiting time distributions u'i(r). If this 
transition is governed by a Poisson process with rate 7 the flow is just this rate 
times the probability that the particle is in state i. Thus, for the transition 
1^2 Ji{t) = 7Pi(t). The probability density P2{t,T) is the flow Ji at time 
t — r times the probability that the process stays in state 2 at least the time 
T. The latter is determined by W2{t). In summary, 



P2{t,T)^^P,{t-T) 



- j dr'w2{r') 




(1) 



By definition we immediately obtain 

P2W = / dr^Piit - r) [1 - y dr'w2{T') 



(2) 



In case of a general waiting time distribution we introduce a waiting time 
dependent conditioned rate ^^(t) from the relation 



T 

Wi{T) = (piir) 1 - y" dr'wiir') 



(3) 



It gives the rate of escaping from state i after having waited the time r and 
determines the flow between i and i + 1 according to Ji{t) — J^dT(f)i{T)Pi{t, r). 
Subsuming we obtain 



OC T 

Pi+i{t) = J drJiit - r) [1 - y dr'wi+iir') 



(4) 



Inserting P2{t,T) from eq.( 1) into this equation for i = 2 finally leads to 



00 00 



Psit) = jdr j dT'-iPi{t - r - t')w2{t') [l - J dT"w^{T" 
00 



(5) 
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The system of cqs. (2) and (5) is completed and becomes closed by the nor- 
malization condition 

P,{t) ^ 1 - P2{t) - Ps{t) . (6) 



Because of the first Poissonian transition the system is generally ergodic and 
possesses one stationary asymptotic solution. 

Further information of the oscillator can be extracted from its spectral proper- 
ties. To this end we consider the stochastic process (output) f{t) as introduced, 
previously. This process corresponds to a sequence of pulses whose lengths p 
are distributed according to Wp(t) = Wsit) and the intervals a between two 
consequent pulses are distributed according to Wa-{t) = /gdr Wi(r)w2(t — r) = 
7 J^dr exp(—'yT)w2{t — r). The spectrum of z{t) — {z{t)) can be calculated 
employing renewal theory [13] and is given by 



^(u) - 4 Re (l-e,M)(l-e.(a;)) 

""^^^ - co^iip) + {a)) 1 - e,(a;)e.(a;) 

where 6p(<^) := (exp(ia;p)) and 0cr(a;) := {exp^icua)) are the characteristic 
functions of the random variables p and a respectively. 

Let us now consider the case where all three transitions are governed by Pois- 
son processes, i.e. W2{t) = 72exp(— 72r) and W3{t) = 73exp(— 73r). After 
taking the derivative of eq. (2) with respect to t we obtain 



^ d 



= - l^rf— 7Pi(i-T)]exp(-72T) = -72P2W+7AW 


A similar calculation can be done for P3{t), leading to 

P3(t) = -73P3(t)+72P2(t), 

i.e. one retrieves the expected rate dynamics. To determine the corresponding 
spectrum as defined in the previous section we use Wp{t) = 73exp(— 73^) and 
'Wait) = :^(exp(-72i) - exp(7i)) which gives (see Fig. (2)) 

^( X 4((7 + 72)^ + ^^) 

^ ' + 4 + :^)((772 + 7273 + 737)^ + u^'h' + 71 + 7I) + ' 

Instead of a fixed rate for the transition 2 — > 3 and 3 — > 1 we now deal with 
the case of a fixed waiting time. The corresponding waiting time distributions 
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are given by W2(t) = S{t — T2) and W3{t) — S{t — T3). Inserting these into 
eqs. (5) and (2) leads to 



T2 T3 

P2{t) = JdTjP,{t-T), Ps{t)= JdTjPi{t-T-T2) 


supplemented by normalization, eq. (6). The calculation of the power spectrum 
is done analogously to the Markovian case, finally giving 



o( N ^ 2cos(a;r3) -2 ^ - 

^"^^ (i + r2 + r3)(a;2 + 72)(i + ^^sin((r2 + r3)a; + 0))' ^""^ 7' 

With T2 fixed to 0.3 and T3 fixed to 0.7 the spectra for different values of the 
rate 7 are presented in figure 2. 
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Fig. 2. The power spectrum for a single stochastic oscillator which Poissonian (left) 
and deterministic waiting time distributions (right). Left: 7 = 0.1,73 = 1. Solid: 
72 = 1, dotted: 72 = 0.1. dashed:72 = 0.01. Right:r2 = 0.3 and n = 0.7. Solid: 
7 = 2.0, dotted:7 = 0.5, dashed: 7 = 0.1. 

The spectra do not differ too much from the Poissonian case. In both cases it 
is possible to have a maximum in the spectrum at frequencies different from 
zero. 



3 Globally coupled three state units 



We now consider a set of N three state oscillators as described in the previous 
section. A global coupling between the single units is introduced by letting the 
rate 7, which governs the transition 1^2, depend on the common output 



N 



(8) 
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where fi{t) is now the output of a single unit. 



In the hmit iV — > cxd the state of the total system can be described by the 
ensemble averaged occupation probabilities Pi{t) that a single oscillator is in 
state i at time t. Since f{t) converges to the occupation probabihty of the 
third state, 7 becomes a function of Psit). 

Equations which govern the dynamics of the ensemble averaged occupation 
probabilities are derived in the same manner as eqs. (5) and (2) except that 
one has to take care that the rate 7 is no longer a constant but is dependent 
on P3{t). One eventually arrives at 



P2{t)^ Jdr-f{Ps{t-r))Pi{t-r)[l- j dr' w^ir')] (9) 


CXD CX) r 

P^it) = JdTj dr' 7(P3(i - T - r'))P^it - r - t')^2(t') [l - / dr" Wsir") 



and normalization, eq.(6). These equations are now nonlinear in the occupa- 
tion probabilities due to the dependency of 7 on the state of the system. 

The system dynamics of the coupled system can be understood by analyzing 
the stationary solutions and its stability properties. The stationary solutions 
P*, P2 and P3* are defined by setting Pi{t) — const, in eq. (9). Using /o°°(l — 
Jgdr' W2/3{T'))dT = J^dr tw2/3{t) =: W2/3, which is the mean waiting time in 
state 2 or 3, respectively, we obtain implicitly 

p: = '-^, pi = ^. (10) 



with T = 7"^(-P3*) +W2+W3 being the mean time of one cycle. In contrast to 
the uncoupled case, where exactly one asymptotic stationary solution exists, 
the asymptotic behavior may now be more complex due to the nonlinearity of 
eqs. (9). Depending on the chosen coupling function there can be several not 
necessarily stable stationary solutions. 

The local stabihty of a stationary solution can be derived by inserting the 
Ansatz Piif) = P* + Ajexp(At), A 7^ 0, into the linearized dynamical equa- 
tions. Therefrom we obtain the characteristic equation 



00 T 

- 1 + 7(P3*) J dre-^^il - J dT'w2{T')) (11) 
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+ (7(^3) -VraA*) JdT[e-'^{l- Jdr'wsir'))] Jdre-'^w^ir) 



In the case with Poissonian transitions the mean waiting times in state 2 and 3 
are given by W2 = I/72 and W3 = I/73. There may be several stable and unsta- 
ble stationary solutions depending on the chosen coupling function. Eq. (11) 
corresponds to the characteristic equation of a two dimensional dynamical 
system: 



+ A(7(P3*) + 72 + 73) + 7(^3*)72 + 7273 + 737(^3*) - 727(^^3 )i'r = 0. 

In order to decide whether there exist limit cycles we calculate the divergence 
of the flow in the phase space which is —{'-/{Ps) +72 + 73). As the rates 7(^3), 
72 and 73 are always positive the divergence is everywhere negative. Therefore 
no limit cycles exist, i.e. the system shows no coherent oscillations. 

For fixed waiting times T2 and T3 the dynamics of the occupation probabilities 
is reduced to 

P^it) = j dT^iPsit - T))P^{t - t) , (12) 


Psit) = J drj{P3{t - T - T2))Pi{t -T-T2), Pi(t) = 1 - P^it) - Ps{t) 


The stationary solutions are given by eq. (10) with W2 = T2 and = T3. The 
characteristic equation is easily derived from the general eq. (11) to be 

7(P3*)(e-^(^+^3) - 1) - 7'(P*)P*(e-^(^+^) - e'^^) - A = 0. (13) 

Such characteristic equations are typical for dynamical systems involving fixed 
delays [14]. 



4 Sigmoidal coupling 



In the following we review the behavior of the coupled system for a sigmoidal 
rate function as presented in [10] 



7(P3(^)) =70 1 + Atanh 



f{t)-r ' 

2(7 



For this rate function, there exists exactly one stationary solution of the system 
dynamics. The coupling strength is mediated by o: The smaller o the stronger 
is the dependence of the rate on Pj,ify- For small a the rate changes rapidly 
from 7i := 7o(l + A) to 72 := 7o(l — A) if the value of P^it) crosses /q. For 
large values of a the rate is nearly independent of P^ii) and equal to 70. 

In case of a Poissonian waiting time distribution this stationary solution is 
necessarily stable as shown in the last section. The high randomness in the 
transitions does not allow the excitation of oscillating behavior, i.e the coupled 
system does not differ substantially from the uncoupled case. 

For fixed waiting times, however, the behavior changes due to the coupling. 
For each value of 70 there is a finite region in which the ensemble oscillates. 
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Fig. 3. Sohitions of the characteristic equation with Re A = for different values of 
7o. The parameter A is fixed to 2/3, /o is chosen to be the stationary solution P| 
at T2 = 0.3, which is about 0.42, 0.16 and 0.04 respectively. The inset shows the 
corresponding frequency at the bifurcation. 

This region grows with increasing values of 70. If cr is chosen too large, coherent 
behavior cannot be observed. Interestingly, there exist values of ti where a 
has to exceed some finite value to trigger the oscillations. We have verified 
this transition between non-oscillating, oscillating and again non-oscillating 
behavior in numerical simulations. 

One may ask, whether there is some connection between the oscillatory be- 
havior of a single oscillator and the oscillatory behavior of the total system. 
To this end we have carried out simulations with the sigmoidal rate function, 
where 70 has been chosen to be 0.06, 0.3 and 1.5 respectively. In the first 
case the rate lies always in the non oscillating region of a single oscillator, in 
the second case, it switches between the oscillating and non oscillating region, 
and in the third case it is always in the oscillating region, a has been fixed 
to 0.001 and the values of /o are chosen to be approximately equal to the 
stationary solution P3*. The results are shown in figure 4. In each case we ob- 
serve coherent oscillations of the total system, although the amplitude of the 
oscillations depends strongly on the chosen rate. Interestingly however, the 
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Fig. 4. Left: Behavior of the ensemble for different values 70 and /o with a 
fixed to 0.001 and A to 2/3 and T2 to 0.3. Solid: 70 = 0.06, /o = 0.04 Dashed: 
7o = 0.3, /o = 0.16 Dotted: 70 = 1-5, /o = 0.42 The simulations have been per- 
formed for a system consisting of 10000 units. Right: The corresponding spectra of 
the common output f{t) 

frequency of the oscillations is nearly independent on 70 although the mean 
time of one cycle of a single oscillator strongly varies with 70. Taking into ac- 
count that the mean waiting time of a Poisson process with rate 7 is 1 /7 the 
mean time for one cycle of a single oscillator lies between 1 + 1/ (70(1 + A)) and 
1 + 1/(70(1 — A)). This corresponds to frequencies between approximately 0.02 
and 0.09 for 70 = 0.06, 0.09 and 0.3 for 70 = 0.3 and 0.3 and 0.7 for 70 = 1.5. 
Obviously the frequencies of the coherent oscillations one observes do not cor- 
respond to these frequencies, especially in the first case they are even of a 
different order of magnitude. However there is a good agreement between the 
observed frequencies and the frequencies at the bifurcation, which are shown 
in the inset of Fig. 3. 



5 Conclusion 

We have investigated stochastic three state oscillators, which we propose as a 
model for excitable systems driven by noise. First we considered the behav- 
ior of the single units and inspected their spectral properties for the case of 
Poissonian and deterministic transitions. It turns out that the process with 
deterministic transitions is a good simple model describing excitable dynamics 
with rich spectral behavior. When coupled globally both models behave quite 
differently. While in the Poissonian case there are no oscillating solutions, in 
the case with fixed waiting times the coupling leads to coherent oscillations of 
the network. A stability analysis reveals, that these oscillations are due to a 
Hopf bifurcation in the dynamics of the occupation probabilities. Interestingly, 
the global oscillations occur even for units which have maximal power at zero 
frequency. It is the common situation of noise induced phase transitions, that 
coupling might induce ordered behavior despite the fact that the constituting 
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elements involve a lot of randomness [2,3,5]. 

This work was supported by DFG-SIId 555. We thank B. Lindner and M. Zaks 
for help and fruitful comments. 
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